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AN EFFICIENT FILTERED SCHEME FOR SOME FIRST ORDER 
HAMILTON-JACOBI-BELLMAN EQUATIONS* 

OLIVIER BOKANOWSKI.t MAURIZIO FALCONE,t SMITA SAHU§ 


Abstract. We introduce a new class of “filtered” schemes for some first order non-linear Hamilton-Jacobi- 
Bellman equations. The work follows recent ideas of Froese and Oberman (SIAM J. Numer. Anal., Vol 51, pp.423- 
444, 2013). The proposed schemes are not monotone but still satisfy some e-monotone property. Convergence results 
and precise error estimates are given, of the order of V Ax where Ax is the mesh size. The framework allows to 
construct finite difference discretizations that are easy to implement, high—order in the domains where the solution 
is smooth, and provably convergent, together with error estimates. Numerical tests on several examples are given 
to validate the approach, also showing how the filtered technique can be applied to stabilize an otherwise unstable 
high-order scheme. 

Key words. Hamilton-Jacobi equation, high-order schemes, e-monotone scheme, viscosity solutions, error 
estimates 

1. Introduction. In this work, our aim is to develop high-order and convergent schemes for 
first order Hamilton-Jacobi (HJ) equations of the following form 

d t v + H{x, Vv) = 0, (t, x) € [0, T\ x R d (1.1) 

n(0, x) = vq(x), x £ R d . (1.2) 


Basic assumptions on the Hamiltonian H and the initial data vq will be introduced in the next 
section. It is well known that, in the one dimensional case, there is a strong link between Hamilton- 
Jacobi equations and scalar conservation laws. Namely, the viscosity solution of the evolutive HJ 
equation is the primitive of the entropy solution of the corresponding hyperbolic conservation law 
with the same hamiltonian. There are several schemes developed for hyperbolic conservation law 
(see references [15] m, 0 , pi). Most of the numerical ideas to solve hyperbolic conservation 
law can be extended to HJ equations. Well known high-order essentially non-oscillatory (ENO) 
scheme have been introduced by A. Harten et al. in m for hyperbolic conservation laws, and then 
extended to HJ equation by Osher and Shu m ■ ENO schemes have shown to have high-order 
accuracy although a precise convergence result is still missing and, for this property, they have 
been quite successful in many applications. Despite the fact that there is no convergence proof of 
ENO schemes towards the viscosity solution of (1.1) in the general case, convergence results may 
hold for related schemes, e.g. MUSCL schemes, as it has been proved by Lions and Souganidis 
in m- Convergence results of some non monotone scheme have also been shown in particular 
cases [5J. Another interesting result has been proved by Fjordholm et al. m, they have shown 
that ENO interpolation is stable but the stability result is not sufficient to conclude total variation 
boundedness (TVB) of the ENO reconstruction procedure. In |T0], a conjecture related to weak 
total variation property for ENO schemes is given. Finally, let us also mention [B], where weighted 
essentially non-oscillatory (WENO) schemes have been applied to HJ equations; the convergence 
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proof of the scheme relies also on the work of Ferretti [9j where higher than first order schemes are 
proposed in a semi-Lagriangian setting, yet with restrictive conditions on the mesh steps. 

In this paper we give a very simple way to construct high-order schemes in a convergent 
framework. It is known (by Godunov’s theorem) that a monotone scheme can be at most first 
order. Therefore it is needed to look for non-monotone schemes. The difficulty is then to combine 


non-monotonicity of the scheme and convergence towards the viscosity solution of (1.1), and also 


to obtain error estimates. In our approach we will adapt a general idea of Froese and Oberman [12] , 
that was presented for stationnary second order Hamilton-Jacobi equations and based on the use 
of a “filter” function. Here we focus mainly on the case of evolutive first order Hamilton-Jacobi 
equation O and an adaptation to the steady case will be also considered. We will design a 
slightly different filter function for which the filtered scheme is still an “e-monotone” scheme (see 
Eq 2. 12|) , but that improves the numerical results. Let us mention also the work [4! for steady 
equations where some e-monotone semi-Lagrangian schemes are studied. 

The paper is organized as follows. In Section 2, we present the schemes and give main conver¬ 
gence results. Section 3 is devoted to the numerical tests on several academic examples to illustrate 
our approach in one and two-dimensional cases. A test on nonlinear steady equations , as well 
an evolutive “obstacle” HJ equation in the form of min(u t + H(x,u x ),u — g(x)) = 0 for a given 
function g are also included in this section. Finally, Section 4 contains our concluding remarks. 

2. Definitions and main results. 

2.1. Setting of the problem. Let us denote by | • | the Euclidean norm on R d (d > 1). The 
following classical assumptions will be considered in the sequel of this paper: 

(Al) vq is Lipschitz continuous function i.e. there exist Lq > 0 such that for every x ,y £ R d , 


(A2) H : 


I Vo (a:) - v 0 {y)\ < L 0 \x - y |. 
satisfies, for some constant C > 0, for all p, q,x,y £ 


( 2 . 1 ) 


I H(y,p) - H(x,p )| < C( 1 + \p\)\y - x\, (2.2) 

and 

I H(x,q) - H(x,p )| < C{ 1 + \x\)\q-p\. (2.3) 


Under assumptions (Al) and (A2) there exists a unique viscosity solution for (1.1) (see Ishii 
mi). Furthermore v is locally Lipschitz continuous on [0,T] x R d . 

For clarity of presentation we focus on the one-dimensional case and consider the following 
simplified problem: 


v t + H(x,v x ) = 0, x £ R, t £ [0,T], (2.4) 

v(0, x) = vo(x), x £ R. (2-5) 

2.2. Construction of the filtered scheme. Let r = At > 0 be a time step (in the form 
of r = ^ for some TV > 1), and Aa; > 0 be a space step. A uniform mesh in time is defined by 
t n := nr , n £ [0, ..., TV], and in space by the nodes Xj := jAx, j £ Z. 

The construction of a filtered scheme needs three ingredients: 

• a monotone scheme, denoted S M 

• a high-order scheme, denoted S A 

• a bounded “filter” function, F : R —> R. 
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The high-order scheme need not be convergent nor stable; the letter A stands for “arbitrary order”, 
following Hlj. For a start, S M will be based on a finite difference scheme. Later on we will also 
propose a definition of S M based on a semi-Lagriangian scheme. 

Then, the filtered scheme S F is defined as 

u n+i = S F (u n )j := S M {u n )j + c T F^ SA ^ i ~ ' SM ^ j ^ , ( 2 . 6 ) 

where e = e T} & x > 0 is a parameter satisfying 

lim e = 0. (2.7) 

( r , Ax )—>-0 

More hints on the choice of e will be given later on. 

The scheme is initialized in the standard way, i.e. 

u°:=vo(xj), Vj e Z. (2.8) 

Now we make precise some requirements on S M , S A and the function F. 

Definition of the monotone finite difference scheme S AI 

Following Crandall and Lions [7{, we consider a finite difference scheme written as u n+1 = S M (u n ) 
with 

S M (u n ){x) := u n {x) - r h M (x, D~u n (x), D+u n {x)), (2.9) 


with 


D±„(x) :=± “ (a ’ ±a A I )~“ (l) , 

where /i M corresponds to a monotone numerical Hamiltonian that will be made precise below. We 
will denote also S M (u n )j := S M (u n )(xj). Therefore the scheme also reads, for all j £ Z, \/n > 0: 


F+ 1 ■_= — t h M { Xj ,D~u^,D + u ”), D ± uf := ± 


u j 

Ax 


( 2 . 10 ) 


(A3) Assumptions on S M 

We will use the following assumptions throughout this paper: 

(i) h M is a Lipschitz continuous function. 

( ii ) (consistency ) Vcc, Vrt, h M (x,u, u) = H(x,u). 

(■ Hi ) (monotonicity ) for any functions it, v, 

u < v => S M {u) < S M {v). 

In pratice condition (A3)-(m) is only required at mesh points and the condition reads 

Uj < Vj, Vj, => S M (u)j < S M (v)j , Vj. (2.11) 

At this stage, we notice that under condition (A3) the filtered scheme is “e-monotone” in the 
sense that 


Uj < Vj , Vj, => S M (u)j < S M (v)j + er, Vj. 
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( 2 . 12 ) 






with e —>• 0 as (r, Ax) — > 0. This implies the convergence of the scheme by Barles-Souganidis 
convergence theorem (see [21 II])- 

Remark 2.1. Under assumption (i), the consistency property ( ii ) is equivalent to say that, 
for any v £ C 2 ([0, T] x R), there exists a constant Cm > 0 independant of Ax such that 


h M (x,D v(x), D + v(x)) — H(x, v x ) 


< C M Ax\\d xx v\ 


(2.13) 


The same statement holds true if (2.13) is replaced by the following consistency error estimate: 
£ S M{v)(t,x) := 


v(t + t,x) — S M (v(t , -))(x) 


- (v t (t,x) + H(x,v x (t,x))) 


< CM\T\\dt t v \\ 00 +Ax\\d xx v\\ 


(2.14) 


Remark 2.2. Assuming ( i), it is easily shown that the monotonicity property (in) is equivalent 
to say that h M = h M (x,u~,u + ) satisfies, a.e. (x,u~ ,u + ) £ R 3 : 


dh 


M 


du~ 


> 0 , 


dh M 

du + 


< 0 , 


and the CFL condition 


t (dh M _ , dh M _ , , 

Ax(du^ {x ^' U ^du^ {x ' U ~' U ) )^ 1 ' 


(2.15) 


(2.16) 


When using finite difference schemes, it is assumed that the CFL condition (2.16) is satisfied, and 
that can be written equivalently in the form 


Co 



where cq is a constant independant of r and Ax. 

Example 2.1. Let us consider the Lax-Friedrichs numerical Hamiltonian is 


•M,LF 


(x, u , u + ) := H(x 


u + w" 1 


\ c o / + 
■)-y(« + - 


(2-17) 


where cq > 0 is a constant. The scheme is consistant; it is furthermore monotone under the 
conditions max XiP \d p H(x,p)\ < cq, and Co^y < 1. 

Definition of the high-order scheme S A : we consider an iterative scheme of “high-order” in the 
form u n+1 = S A (u n ), written as 

S A (u n )(x) = u n (x ) — rh A (x, D k ’~u n (x),..., D~u n (x), D + u"(x),..., D k ' + u n (x)), 


where h A corresponds to a “high-order” numerical Hamiltonian, and H f,± u( x) := 1 

for £ = 1,..., k. To simplify the notations we may write (|2. 181 in the more compact form 


S A (u n )(x) = u n (x) - rh A (x, D ± u r 


0 * 0 ) 


(2.18) 


even if there is a dependency on £ in (D l ^u n (x))z = 

(A4) Assumptions on S A : 
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We will use the following assumptions: 

(i) h A is a Lipschitz continuous function. 

(ii) (high-order consistency) There exists k > 2, for all £ £ [1,. ..,k], for any function v = 
v(t,x) of class C l+1 , there exists Ca,i > 0, 


£s A (v)(t,x) := 


v(t + r, x) - S A (v(t , .))(x) 


- (v t (t,x) + H(x,v x (t,x))) 


< C^fr'liaf+^IU+A^II^II 


(2.19) 

( 2 . 20 ) 


Here v^ denotes the £-th derivative of v w.r.t. x. 

Remark 2.3. The high-order consistency implies, for all £ £ [1, ..., k], and for v £ C f+1 (R), 


h A (x,...,D v, D + v ,...) -H(x,v x ) 


< CaA^vWoo^. 


Example 2.2. (Centered scheme) A typical example with k = 2 is obtained with the centered 
TVD (Total Variation Diminishing) approximation in space and the Runge-Kutta 2nd order scheme 
in time (or Heun scheme): 


u n — u n 

S 0 {u n ) 3 := u] - rH{x jt J+1 2Ax J ~ 1 ), (2.21a) 

and 

S A (u):=^(u + S 0 (S 0 (u))). (2.21b) 

Of course there is no reason for the centered scheme to be stable (as it will be shown in the 
numerical section). Using a filter will help stabilize the scheme. A similar example with k = 3 
can be obtained with any third order finite difference approximation in space and the TVD-RK3 
scheme in time un- 

Definition of the filter function F. We recall that Froese and Oberman’s filter function used in |12| 
is: 


F(x) = sign(x) max(l — ||x| — 11,0) 


x |x| < 1. 

0 jxj > 2 . 

—x + 2 1 < x < 2. 

—x — 2 —2 < x < —1. 


In the present work we define a new filter function simply as follows: 


F(x) := x1 n <;l 


x if |x| < 1, 

0 otherwise. 


( 2 . 22 ) 


The idea of the present filter function is to keep the high-order scheme when \h A — h M \ < e 
(because then \S A — S M \/(re) < 1 and S F = S M + reF( sA f^ M ) = S A ), whereas F = 0 and 
S F = S 1 '' 1 if that bound is not satisfied, i.e., the scheme is simply given by the monotone scheme 
itself. Clearly the main difference is the discontinuity at x = — 1,1. 
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Fig. 2.1. Froese and Oberman's filter (left), new filter (right) 


2.3. Convergence result. The following theorem gives several basic convergence results for 
the filtered scheme. Note that the high-order assumption (A4) will not be necessary to get the 
error estimates ( [i)-(ii )■ It will be only used to get a high-order consistency error estimate in the 
regular case (part (in)). Globally the scheme will have just an 0(VAx) rate of convergence for 
just Lipschitz continuous solutions because the jumps in the gradient prevent high-order accuracy 
on the kinks. 

Theorem 2.1. Assume (A1)-(A2), and let vq be bounded. We assume also that S M satisfies 
(A3), and |F| < 1. Let u n denote the filtered scheme (2.6). Let v ” := v(t n , Xj) where v is the exact 
solution of (2.4). Assume 


0 < e < cq VAx 


(2.23) 


for some constant cq > 0. 

(i) The scheme u n satisfies the Crandall-Lions estimate 


|u"-V B ||oo 


< C\/Ax, V n = 0, ...,N. 


(2.24) 


for some constant C independent of Ax. 

(ii) (First order convergence for classical solutions.) If furthermore the exact solution v belongs 
to C 2 ([0,T] x M), and e < cqAx (instead of ( 2.23\ ), then, we have 


\u n -v n \\ 00 <CAx, n = 0,..., N, 


(2.25) 


for some constant C independent of Ax. 

(in) (Local high-order consistency.) Assume that S A is a high-order scheme satisfying (A)) for 
some k > 2. Let 1 < I < k and v be a C e+1 function in a neighborhood of a point (t, x) £ (0, T) x M. 
Assume that 


(2.26) 


(Cj 4,1 + C M ) ^||«tt||oo T + ll^xlU A.t ) < e. 
Then, for sufficiently small t n — t, Xj — x, t, Ax, it holds 


S F (v n )i = S A (v n ) 3 
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and, in particular, a local high-order consistency error for the filtered scheme S F holds: 


£ SF (v n ) j = E SA (v n ) j = 0(Ax e ) 


(the consistency error £$ A defined in (2.19)/ 

Proof, (i) Let w™ +1 = S M (w n )j be defined with the monotone scheme only, with 
Vo(xj) = vPy By definitions, 

u] +1 - w” +1 = S M (u n )j - S M (w n )j + ctF(.) 

Hence, by using the monotonicity of S M , 

max |it" +1 — wf +1 1 < max |it" — w" | + er, 

3 J 

and by recursion, for n < N, 


max | uf — wf \ < eur < Te. 

3 


On the other hand, by Crandall and Lions [7|, an error estimate holds for the monotone scheme: 


max | w™ — u”| < CVAx, 

3 

for some C > 0. By summing up the previous bounds, we deduce 

max |< - vf \ < CVA^c + Te, 

3 

and together with the assumption on e, it gives the desired result. 

v n+l _ S M (v n ) 0 

( ii ) Let £” := —-. If the solution is C l regular with bounded second order 

J r 

derivatives, then the consistency error is bounded by 


\£?\<C m {t + Ax). 

Hence 


(2.27) 


l^n+l _ v n+l | | o M/„.n\ qM 


= \S M (u n ) j -S M (v n ) j + T£? + reF(.)\ 


< IK-^ n IU + r|r 


By recursion, for nr < T , 


\\u n - «"||oo < ||w° - u°||<x> +T( max ||£ ||oo + e)- 

0<fe<JV-l 


Finally by using the assumption on e, the bound (2.27) and the fact that r = 0(Ax) (using CFL 
condition (2.17)), we get the desired result. 

(Hi) To prove that S F (v n )j = S A (v n )j, one has to check that 


|£>"),--S Af (t;") J -| 


< 1 
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as (t, Ax) —>• 0. By using the consistency error definitions, 


| S A (v n )j-S M (y n ) 


n+1 




T 

n+1 


+ v t (t n ,Xj) + H(xj,v x (t n ,Xj)) 


- S A {v n ) j 


+ v t (t n ,Xj) + H(xj,v x (t n ,Xj)) 


< \£ s u(v n )j \ + \£ S A(v n )j\ 

< (C A , i + Cm){j\\vu\\oo + Ax||u x 


0 


Hence the desired result follows. □ 

Remark 2.4. Other approaches. It is already known from the original work of Osher and 
Shu that it is possible to modify an ENO scheme in order to obtain a convergent scheme. For 
instance, if D ± A uf denotes a high-order finite difference derivative estimate (of ENO type), a 
projection on the first order finite difference derivative D ± uf can be used, up to a controlled error 
(see in particular Remark 2.2 of [2Vf): 


instead of D ± ' A uf, 

where P[ a ,b\{y ) the projection defined by: 


use ,MAx](P ’ Uj ) 


if a — b<y<a + b 


P[a,b] (y) ■= \ a — b if y < a- b 
a + b if y > a + b 

and M > 0 is some constant greater than the expected value \\u xx (t n ,Xj)\. However, we emphasize 
that in our approach we do not consider a projection but a perturbation with a filter, which is sligthly 
different. Indeed, by using a projection into an interval of the form [a — MAx, a + M Ax] where 
a = D ± u™, numerical tests show that we may choose too often one of the extremal values a±MAx 
which is then produces an overall too big error (worse than using the first order finite differences). 
Following the present approach, we would rather advice to use, 

instead of D ± ' A uf, the value P[D±u n ,MAx](D ± ' Au 'j) 
where P[ a ,b]{y ) defined by: 


P[a,b]{y) ■= 


y if a — b < y < a + b 

a tf y [a ~ b, a + b] 


Remark 2.5. Filtered semi-Lagrangian scheme. Let us consider the case of 

H(x,p) := minmax{—/(x, a, b).p — £(x, a, &)}, (2.28) 

b£B cl€lA 

where A c R m and B C R n are non-empty compact sets (with m,n > 1), f : x A x B —> 

and l : M. d x A x B —► R. are Lipschitz continuous w.r.t. x: 3 L > 0, V(a, b) £ A x B, Vx, y: 

max(|/(x, a, b) — f(y, a, b)\, \I(x, a, b) — I(y, a, 6)|) < L\x - y\. (2.29) 

(We notice that (A2) is satisfied for hamiltonian functions such as flA28l ).j Let [w] denote the 
P 1 -interpolation of u in dimension one on the mesh (Xj), i.e. 


x € \Xj,X j+ 1 ] 


r i/ \ x j +1 x 

[“](*) : = Ax U 3 


X — Xj 

Ax 


u j+ 1- 


(2.30) 











Then a monotone SL scheme can be defined as follows: 


S M (u n )j := min max ( [u n ](xj + rf(xj,a,b)) + rl(xj, a, b) ). (2-31) 

adA b£B \ ' J 


A filtered scheme based on SL can then be defined by rt2.6k and (2.31). Convergence result as well 


as error estimates could also be obtained in this framework. (For error estimates for the monotone 
SL scheme, we refer to mw 


2.4. Adding a limiter. The basic filter scheme (2.6) is designed to be of high-order where 


the solution is regular and when there is no viscosity aspects. However, for instance in the case 
of front propagation, it can be observed that the filter scheme may let small errors occur near 
extrema, when two possible directions of propagation occur in the same cell. 

This is the case for instance near a minima for an eikonal equation. In order to improve the 
scheme near extrema, we propose to introduce a limiter before doing the filtering process. Limiting 
correction will be needed only when there is some viscosity aspect (it is not needed for advection). 


Let us consider the case of front propagation, i.e., equation of type (2.4), now with 


H(x, v x ) = max (f(x, a) 

aGA 


(2.32) 


(i.e., no distributive cost in the Hamiltonian function). 

In the one-dimensional case, a viscosity aspect may occur at a minima detected at mesh point 
Xi if 


min a f(xj, a) < 0 and max a f(xj, a) > 0. (2.33) 

In that case, the solution should not go below the local minima around this point, i.e., we want 

u] +1 > u min j := min(u"_ x , it", Uj +1 ), (2.34) 

and, in the same way, we want to impose that 

u" +1 < u max j := max(u"_ 1; u”, u" +1 ). (2.35) 

If we consider the high-order scheme to be of the form u" +1 = uf — rh A (u n ), then the limiting 
process amounts to saying that 


h A (u n )j < h r j 


and 


. u- — u mnr A 

h A (u n )j > h™ m := -A - CCfLA, 


This amounts to define a limited h A such that, if (2.33) holds at mesh point x,, then 


h A (u n )j := min ( ma x(h A (u n )j, h r f lm ), hL 


and, otherwise, 
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Then the filtering process is the same, using h A instead of h A for the definition of the high-order 
scheme S F . 

For two dimensional equations a similar limiter could be developped in order to make the 
scheme more efficient at singular regions. However, for the numerical tests of the next section (in 


two dimensions) we will simply limit the scheme by using an equivalent of (2.34)-(2.35). Hence, 
instead of the scheme value u™ +1 = S A (u n )ij for the high-order scheme, we will update the value 
by 


u™ + = min(max(5‘ (u n )ij , uf m ), u 


max\ 


(2.36) 


where vf in = min(u£, < ±1J , u” j±1 ) and uf ax = max(u£, u? ±1J , K J±1 ). 


2.5. How to choose the parameter e: a simplified approach. The scheme should switch 
to high-order scheme when some regularity of the data is detected, and in that case we should 
have 


S A (v)-S M (v) 


h A (.)-h M (-) 

er 


e 


In a region where a function v = v(x) is regular enough, by using Taylor expansions, zero order 
terms in h A (x, D±v) and h M {x,D ± v) vanish (they are both equal to H(x,v x (x))) and it remains 
an estimate of order Ax. More precisely, by using the high-order property (A4) we have 

h A (xj, D ± Vj) = H(xj,v x (xj)) + 0{Ax 2 ). 


On the other hand, by using Taylor expansions, 

Dvf = v x (xj) ± -v xx (xj)Ax + 0(Ax 2 ), 

Hence, denoting h M = h M (x,u~,u + ), it holds at points where h M is regular, 


h M (. xj , Dv, , DvJ ) = H(xj , v x (. 


O')) + 2 Vx 


■(Xj) 


dhf 




\ du + du~ 


■ 


+ 0( Ax 2 ). 


Therefore, 


\h A (v)-h M (v)\ 


2 \ v xx{ x j)\ 


dhf 

lhj+ 


dhf 

du~ 


Ax + 0( Ax 2 ). 


Hence we will make the choice to take e roughly such that 


1 

2 


\ v xx( x j)\ 


dhf 

du+ 


dhf 

du~ 


Ax < e 


(2.37) 


(where hf = h M (xj, v x (xj), v x (xj))). Therefore, if at some point Xj (2.37) holds, then the scheme 
will switch to the high-order scheme. Otherwise, when the expectations from h M and h A are 
different enough, the scheme will switch to the monotone scheme. 

In conclusion we have upper and lower bound for the switching parameter e: 

• Choose e < CqVAx for some constant Co > 0 in order that the convergence and error 
estimate result holds (see Theorem 2.1). 
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• Choose e > CiAx, where C\ is sufficiently large. This constant should be choosen roughly 
such that 


1 


dh M . . dh M . 

- g^(.,v x ,v x ) 


< Cl. 


where the range of values of v x and v xx can be estimated, in general, from the values of 
(t’o)xx and the Hamiltonian function H. Then the scheme is expected to switch to 
the high-order scheme where the solution is regular. 

3. Numerical tests. In this section we present several numerical tests in one and two dimen¬ 
sions. Unless otherwise indicated, the filtered scheme will refer to the scheme where the high-order 


Hamiltonian is the centered scheme in space (see Remark 2.2), with Heun (RK2) scheme discretiza¬ 
tion in time (see in particular Eqs. (2.21a)-(2.21b)). Hereafter this scheme will be referred as the 
“centered scheme”. 

The monotone finite difference scheme and function h M will be made precise for each example. 

For the filtered scheme, unless otherwise precised, the switching coefficient e = 5Ax. will be 
used. In practice we have numerically observed that taking e = CiAx with c\ sufficiently large does 
not much change the numerical results in the following tests. All the tested filtered schemes (apart 
from the steady and obstacle equations) enters in the convergence framework of the previous 
section, so in particular there is a theoretical convergence of order \JAx under the usual CFL 
condition. 

In the tests, the filtered scheme will be in general compared to a second order ENO scheme 
(for precise definition, see Appendix[A|, as well as the centered (a priori unstable) scheme without 
filtering. 

In several cases, local error in the L 2 norms are computed in some subdomain D, which, at a 
given time t n , corresponds to 


e,L? := 


E i 

{z, XiGD} 


'(tn : %i) 


1/2 


The first two numerical examples deal with one-dimensional HJ equations, examples 3 and 
4 are concerned with two-dimensional HJ equations, and the last three examples will concern a 
one-dimensional steady equation and two nonlinear one-dimensional obstacle problems. 


Example 1. Eikonal equation. We consider the case of 

Ut + Kl=0, t £ (0, T), x £ (—2, 2), 

u(0, x) = vq(x) := max(0,1 — x 2 ) 4 , x £ (—2, 2). 


(3.1) 

(3.2) 


In Table 3.1 we compare the filtered scheme (with e = 5Ax) with the centered scheme and the 


ENO second order scheme, with CFL = 0.37 and terminal time T = 0.3. For the filtered scheme, 
the monotone hamiltonian used is h M (x, v~, v + ) := max(v~, —v + ). 

As expected, we observe that the centered scheme alone is unstable. On the other hand, the 
filtered and ENO schemes are numerically comparable in that case, and second order convergent 
(the results are similar for the L 1 and the L°° errors). 

Then, in Table pT~2| we consider the same PDE but with the following reversed initial data: 


Vq( x ) := ~ max(0,1 — x“) , x £ (—2, 2). 


(3.3) 


In that case the centered scheme alone is unbounded. The filtered scheme (with e = 5Ax) is second 


order. However, here, the limiter correction as described in section (2.4) was needed in order to 
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filtered (e 

= 5Az) 

centered 

EN02 

M 

N 

L 1 error 

order 

L 2 error 

order 

L 2 error 

order 

40 

9 

7.51E-03 

- 

1.18E-01 

- 

1.64E-02 

- 

80 

17 

3.36E-03 

1.16 

1.14E-01 

0.06 

4.38E-03 

1.91 

160 

33 

8.02E-04 

2.07 

1.13E-01 

0.00 

1.19E-03 

1.87 

320 

65 

1.80E-04 

2.16 

1.13E-01 

0.00 

3.22E-04 

1.89 

640 

130 

4.53E-05 

1.99 

1.13E-01 

0.00 

8.22E-05 

1.97 


Table 3.1 

(Example 1 with initial data {3.2) ) L 2 errors for filtered scheme, centered scheme, and ENO second order 
scheme 



filtered (e 

= 5Ax) 

centered 

EN02 

M 

N 

error 

order 

error 

order 

error 

order 

40 

9 

1.27E-02 

- 

2.03E-02 

- 

2.60E-02 

- 

80 

17 

3.17E-03 

2.00 

8.96E-03 

1.18 

8.00E-03 

1.70 

160 

33 

7.90E-04 

2.01 

1.06E-02 

-0.24 

2.50E-03 

1.68 

320 

65 

1.97E-04 

2.00 

1.26E-01 

-3.57 

7.80E-04 

1.68 

640 

130 

4.92E-05 

2.00 

1.06E+02 

-9.71 

2.44E-04 

1.67 


Table 3.2 

(Example 1 with initial data ( |3.3[ ) L 2 errors for filtered scheme, centered scheme, and ENO second order 

scheme. 


get second order behavior. We also observe that the filtered scheme gives better results than the 
ENO scheme. (We have also numerically tested the ENO scheme with the same limiter correction 
but it does not improve the behavior of the ENO scheme alone). 

In conclusion, this first example shows firstly, that the filtered scheme can stabilize an otherwise 
unstable scheme, and secondly that it can give the desired second order behavior. 



Fig. 3.1. (Example 1) With initial data ( |3.2| ) (left), and plots at time T = 0.3 with centered scheme - middle 
- and filtered scheme - right, using M = 160 mesh points. 


Example 2. Burger’s equation. 

In this example an HJ equivalent of the nonlinear Burger’s equation is considered: 

v t + ] ) \v x \ 2 =0, t> 0, (-2,2) 

u(0, x) = vq{x) := max(0,1 — x 2 ), x £ (—2,2) 
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(3.4a) 

(3.4b) 




































Fig. 3.2. (Example 1) With initial data ( |3.3[ ) (left), and plots at time T = 0.3 with centered scheme - middle 
- and filtered scheme - right, using M = 160 mesh points. 



filtered (e 

= 5Az) 

centered 

EN02 

M 

N 

error 

order 

error 

order 

error 

order 

40 

9 

2.06E-02 

- 

2.07E-02 

- 

2.55E-02 

- 

80 

17 

6.24E-03 

1.73 

6.24E-03 

1.73 

8.24E-03 

1.63 

160 

33 

1.85E-03 

1.76 

1.85E-03 

1.76 

2.81E-03 

1.55 

320 

65 

5.51E-04 

1.74 

5.51E-04 

1.74 

1.03E-03 

1.45 

640 

130 

1.68E-04 

1.71 

1.68E-04 

1.71 

3.74E-04 

1.47 


Table 3.3 

(Example 2) L 2 errors for filtered scheme, centered scheme, and ENO second order scheme. 


with Dirichlet boundary condition on (—2,2). Exact solution is knownQ In order to test high 
order convergence we have considered the smoother initial data which is the one obtained from 
(3.41 at time to '■= 0.1, i.e. : 


w t + )fw x \ 2 = 0, t > 0, x G (-2,2). 
w(0,x) := v(to,x), x £ (—2,2), 


(3.5a) 

(3.5b) 


with exact solution w(t, x ) = v(t + to, x). 

An illustration is given in Fig. |3.3[ For the filtered scheme, the monotone hamiltonian used is 
h M (x,v~,v+) := ±( v~) 2 l„- >0 + lu+<o- 


Errors are given in Table (3.3), using CFL=0.37 and terminal time T = 0.3. 


In conclusion we observe numerically that the filtered scheme keeps the good behavior of the 
centered scheme (here stable and almost second order). 

Example 3. 2D rotation. We now apply filtered scheme to an advection equation in two 


It holds 


v(t, x) = 


(max(0,1 — |#|))^ 
2 1 


!{*>!} + 


(1 - 2 t ) 2 - | x | 2 
1 - 2 1 
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Fig. 3.3. (Example 2) Plots at t = 0 and t = 0.3 with the filtered scheme. 


dimensions: 


v t - yv x + xv y = 0, (x,y) £ fl, t > 0, (3.6) 

I _ ^ _ I') 2 _ y 2 

u(0, x, y) = Vq(x, y ) := 0.5 — 0.5 max(0,-^-) 4 (3-7) 

1 ~ r o 

where fl := (—A, A) 2 (with A = 2.5), ro = 0.5 and with Dirichlet boundary condition v(t, x) = 0.5, 
x £ dQ,. This initial condition is regular and such that the level set Vo(x,y) = 0 corresponds to a 
circle centered at (1,0) and of radius r 0 . 

In this example the monotone numerical Hamiltonian is defined by 

h M (u~,u+,u~,u+) := max(0, /i(a, x, y))u~ + min(0, f t (a, x, y))u+ (3.8) 

+ max(0 ,f 2 (a,x,y))u~ + min(0, f 2 (a, x, y))u+ 


and the high-order scheme is the centered finite difference scheme in both spacial variables, and 


RK2 in time. The filtered scheme is otherwise the same as (2.6). However it is necessary to use a 
greater constant C\ is the choice e = ciAx in order to get (global) second order errors. We have 
used here e = 20Aa;. 

On the other hand the CFL condition is 


^ * = C ° ( Aa: + ' 


(3.9) 


where here Cq = 2.5 (an upper bound for the dynamics in the considered domain Q). In this test 
the CFL number is fi := 0.37. 


Results are shown in Table 3.4 for terminal time time T := tt/2. Although the centered scheme 


is a priori unstable, in this example it is numerically stable and of second order. So we observe 
that the filtered scheme keep this good behavior and is also or second order (ENO scheme gives 
comparable results here). 

Example 4. Eikonal equation. In this example we consider the eikonal equation 


v t + |Vu| = 0, (a :,y) £ H, t > 0 
in the domain H := (—3, 3) 2 . The initial data is defined by 
vo (x,y) = 


(3.10) 

(3.11) 


0.5 — 0.5 max ( max(0 


1 - (x - l) 2 - y 2 
i-r 2 0 


max(0 


i- (s + i f-y 2 

l-rl 
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filtered 

centered 

ENO 

Mx 

Ny 

L 2 error 

order 

L 2 error 

order 

L 1 error 

order 

20 

20 

5.05E-01 

- 

5.05E-01 

- 

6.99E-01 

- 

40 

40 

1.48E-01 

1.77 

1.48E-01 

1.77 

4.66E-01 

0.58 

80 

80 

3.77E-02 

1.98 

3.77E-02 

1.98 

2.04E-01 

1.19 

160 

160 

9.40E-03 

2.00 

9.40E-03 

2.00 

5.50E-02 

1.89 

320 

320 

2.34E-03 

2.01 

2.34E-03 

2.01 

1.29E-02 

2.10 


Table 3.4 

(Example 3) Global L 2 errors for the filtered scheme, centered and second order ENO schemes (with CFL 0.37 )■ 



Fig. 3.4. (Example 3) Filtered scheme, plots at time t = 0 (left) and t = 7t/2 (rigth) with M = 80 mesh points. 


The zero-level set of vq corresponds to two separates circles or radius r 0 and centered in A = (1,0) 
and B = (—1,0) respectively. Dirchlet boundary conditions are used as the previous example. 
The monotone hamiltonian h M used in the filtered scheme is in Lax-Friedrichs form: 




h (x,tq ,uj,u 2 ,uj) = H(x, 


~) 


c. 


c„ 


-7T ( U 1 “ u l ) - ( U 2 - U 2 )> 


(3.12) 


where, here, 


C X =Cy = 


1. We used the CFL condition /.i = 0.37 as in (3.9). Also, the simple 


limiter (2.36) was used for the filtered scheme as described in Section 2.4 It is needed in order to 


get a good second order behavior at extrema of the numerical solution. The filter coefficient is set 
to e = 20Ax as in the previous example. 

Numerical results are given in Table 3.5 showing the global L 2 errors for the filtered scheme, 
the centered scheme, and a second order ENO scheme, at time t = 0.6. We observe that the 
centered scheme has some unstabilities for fine mesh, while the filtered performs as expected. 

Example 5 Steady eikonal equation. We consider a steady eikonal equation with Dirichlet 
boundary condition, which is taken from Abgrall pQ: 


\vx\ = f(x) x £ (0,1), 
v(0) = v(l) = 0, 

y^2+2 


(3.13a) 

(3.13b) 


where /( x) = 2>x 2 + a, with a = 2 x 0 -i an< ^ x ° = 4 ^T' -^ xac ^ solution is known: 

v(x) := 


x + ax 

1 + a — ax — x 3 


x € [0, a? 0 ], 
x G [x 0 ,1]. 
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(3.14) 






















































filtered (e 

= 20Ax) 

centered 

EN02 

Mx 

Ny 

L 1 error 

order 

L 1 error 

order 

L 2 error 

order 

25 

25 

5.39E-01 

- 

6.00E-01 

- 

5.84E-01 

- 

50 

50 

1.82E-01 

1.57 

2.25E-01 

1.41 

2.11E-01 

1.47 

100 

100 

3.72E-02 

2.29 

8.46E-02 

1.41 

6.88E-02 

1.62 

200 

200 

9.36E-03 

1.99 

3.53E-02 

1.26 

2.02E-02 

1.76 

400 

400 

2.36E-03 

1.99 

1.36E-01 

-1.95 

5.98E-03 

1.76 


Table 3.5 

(Example 4) Global L 2 errors for filtered scheme, centered and second order ENO schemes. 



Fig. 3.5. (Example 4) Plots at times t = 0 (top) and t = tt/ 2 (bottom) for the filtered scheme with M = 50 
mesh points. The figures to the right represent the 0-level sets. 


The steady solution for (3.13) can be considered as the limit lim v(t,x) where v is the solution of 

- t—> oo 

the time marching corresponding form: 

v t + \v x \ = f(x) x G (0,1), t > 0, 
v(t, 0) = v(t, 1) = 0, t > 0. 

In this example the upwind monotone scheme is used: 


(3.15a) 

(3.15b) 


u n+1 - u n 

h M {) j := J J 


- max { Ax . aJ } - 


Ax 


the high-order scheme will be the centered scheme, and the filtered scheme ( |2.6| ) will be used with 
e = 5Ax. The iterations are stopped when the difference between too successive time step is small 
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enough or a fixed number of iterations is passed, i.e., in this example, 

||u" +1 — u n \\i,°° '■= max|u” +1 — <| < 10 -6 or n > N max := 5000. (3.16) 

i 

As analyzed in [4j for e-monotone schemes, for a given mesh step, even if the iterations may not 
converge (because of the non monotony of the scheme), it can be shown to be close to a fixed point 
after enough iterations. 



filtered 

centered 

filtered + ENO 

M 

error 

order 

error order 

error 

order 

50 

2.16E-03 

- 

NaN 

5.29E-03 

- 

100 

7.14E-04 

1.60 

NaN 

1.35E-03 

1.97 

200 

2.17E-04 

1.72 

NaN 

3.42E-04 

1.98 

400 

6.32E-05 

1.78 

NaN 

8.61E-05 

1.99 

800 

2.17E-05 

1.54 

NaN 

2.16E-05 

2.00 


Table 3.6 


(Example 5) Global errors for filtered scheme, compared with the centered (unstable) scheme, and a filtered 
ENO scheme. 



Example 6 Advection with an obstacle. Here we consider an obstacle problem, which is taken 
from [3]: 


min(vt + v x , v — g(x )) =0, t > 0, x € [—1,1], (3-17) 

Fig. 3.6. ( Exarrljilic 1 )) Eifie&ect cfs \: 1/1 \t ion, with M = 50 mesh points. (3.18) 

together with periodic boundary condition. The obstacle function is g(x) := sinlfKx). In this case 
exact solution is given by: 

! max(uo(a: — at), g(x)) if t < i 

max(t) 0 (a: - at), g(x), -l xe [ 0 . 5> i]) if £ e [i, §], (3.19) 

max(co (x-at),g{x),l x€ [ _ 1;t _| ]u[0 . 5 5 i ] ) if t G [|,1], 

Results are given in Table [3771 for terminal time t = 0.5. Errors are computed away from singular 
points, i.e., in the region [—1,1] \ ( Uj = i i3 [s* — <5, Sj + <5]) (where si = —0.1349733, s 2 = 0.5 and 
s 3 = 2/3 are the three singular points. Filtered scheme is numerically of second order (ENO gives 
comparable results here). 
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Errors 

filtered e 

= 5Az 

centered 

EN02 

M 

N 

error 

order 

error 

order 

error 

order 

40 

20 

7.93E-03 

2.03 

1.63E-02 

1.54 

2.14E-02 

1.59 

80 

40 

1.84E-03 

2.10 

2.98E-02 

-0.87 

7.75E-03 

1.46 

160 

80 

3.92E-04 

2.24 

1.46E-02 

1.03 

1.07E-03 

2.86 

320 

160 

9.67E-05 

2.02 

8.02E-03 

0.86 

2.72E-04 

1.97 

640 

320 

2.40E-05 

2.01 

4.10E-03 

0.97 

6.92E-05 

1.98 


Table 3.7 

(Example 6) L°° errors away from singular points, for filtered scheme, centered scheme, and second order 
ENO scheme. 



Fig. 3.7. (Example 6) Plots at T=0(initial data), T=0.3, T=0.5. 


Example 7 Eikonal with an obstacle. We consider an Eikonal equation with an obstacle term, 
also taken from [3]: 


min(v t + \v x \, v — g(x)) = 0, t > 0, x £ [— 1,1], (3.20) 

uo(a;) = 0.5 + sinfax) a; £ [—1,1], (3-21) 

with periodic boundary condition on (—1,1) and g{x) = sin(irx). In this case the exact solution 
is given by: 


where v is the solution of the 
holds, which simplifies to 

v(t, x) := < 


v(t, x) = ma x(v(t, 

x),g(x)). 

(3.22) 

Eikonal equation v t +\v x \ 

= 0. The formula v(t, x) = min yg r,,,_ 

- t,x-\-t] 

v 0 (x + 1) 

if x < —0.5 — t 


-0.5 

if a; £ [—0.5 — t, —0.5 + t], 

(3.23) 

min(u 0 (a; — t), v 0 (x + t)) 

if x > —0.5 + t, 



Results are given in Table [3~8] for terminal time T = 0.2. Plots are also shown in Figure [3i~8] 
for different times (for t > ^ solution remains unchanged). 

4. Conclusion. We propose a “filtered” scheme which behaves as a high-order scheme when 
the solution is smooth and as a low order monotone scheme otherwise. It has a simple presentation 
that is easy to implement. Rigorous error bounds hold, of the same order as the Crandall-Lions 
estimates in vAx where Ax is the mesh size. In the case the solution is smooth a high-order 
consistency error estimate also holds. We show on several numerical examples the ability of the 
scheme to stabilize an otherwise unstable scheme, and also we observe a precision similar to a 
second order ENO scheme on basic linear and non linear examples. 
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Errors 

filtered 

EN02 

M 

error 

order 

error 

order 

40 

3.74E-03 

- 

6.85E-03 

- 

80 

6.26E-04 

2.58 

2.12E-03 

1.69 

160 

1.13E-04 

2.47 

6.80E-04 

1.64 

320 

2.26E-05 

2.32 

2.18E-04 

1.64 

640 

5.50E-06 

2.04 

6.96E-05 

1.65 


Table 3.8 

Filtered scheme and ENO scheme at time t = 0.2 



Fig. 3.8. (Example 1) Plots at times t = 0, t = 0.2 and t = 0.4. The dark line is the numerical solution, 
similar to the exact solution, and the ligth line is the obstacle function. 


On going works concern the application of the present approach to some front propagation 
equations. 

Appendix A. An essentially non-oscillatory (ENO) scheme of second order. 

We recall here a simple ENO method of order two based on the work of Osher and Shu 2Q\ 

for Hamilton Jacobi equation (the ENO method was designed by Harten et al. lT7j for the approx¬ 

imation solution of non-linear conservation law). 

Let m be the rninmod function defined by 

f o if |a| < |6|, ab > 0 

m(a, b) = < b if |6| < |6|, ab > 0 (A.l) 

[ 0 if ab < 0 

(other functions can be considered such as m{a , b) = a if |a| < |6| and m(a, b) = b otherwise). Let 
D ± Uj = ±(uj±i — Uj)/Ax and 

jj 2 U j +1 ~ ^Uj + Uj- 1 

J As 2 

Then the right and left ENO approximation of the derivative can be defined by 

D ± Uj = D ± Uj =F -Air m{D 2 Uj,D 2 Uj±\) 

and the ENO (Euler forward) scheme by 

So(u)j := Uj — rh M (xj, D~Uj, D + Uj). 

The corresponding RK2 scheme can then be defined by S(u) = \{u + Sq(Sq(u))). 
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